Tugas Akhir TPM

Aida Darajati

16 May, 2024

1 Data

1.1 Input Data

# Ganti folder
path <- function() gsub ("\\\\", "/", readClipboard() )
setwd("C:/Users/Fathan/Documents/Obsidian Vault/2. Kuliah/Smt 6/Tugas TPM Aida")

# Baca Data
library(dplyr)
raw.data <- read.csv("Data.csv", sep=";") %>% #Baca csv
  as.data.frame() # Ubah sebagai dataframe
data <- raw.data

# Tampilkan data
library(DT)
datatable(data)

1.2 Check Data

Cek tipe data

str(data)

## 'data.frame':    7385 obs. of  11 variables:
##  $ CO2.Emissions.g.km.             : int  196 221 136 255 244 230 232 255 267 212 ...
##  $ Engine.Size.L.                  : int  2 2 2 4 4 4 4 4 4 2 ...
##  $ Fuel.Consumption.City..L.100.km.: int  10 11 6 13 12 12 12 13 13 11 ...
##  $ Fuel.Consumption.Hwy..L.100.km. : int  7 8 6 9 9 8 8 9 10 8 ...
##  $ Fuel.Consumption.Comb..L.100.km.: int  9 10 6 11 11 10 10 11 12 9 ...
##  $ Make                            : chr  "ACURA" "ACURA" "ACURA" "ACURA" ...
##  $ Model                           : chr  "ILX" "ILX" "ILX HYBRID" "MDX 4WD" ...
##  $ Vehicle.Class                   : chr  "COMPACT" "COMPACT" "COMPACT" "SUV - SMALL" ...
##  $ Cylinders                       : int  4 4 4 6 6 6 6 6 6 4 ...
##  $ Transmission                    : chr  "AS5" "M6" "AV7" "AS6" ...
##  $ Fuel.Type                       : chr  "Z" "Z" "Z" "Z" ...

Cek data hilang

colSums(is.na(data))

##              CO2.Emissions.g.km.                   Engine.Size.L. 
##                                0                                0 
## Fuel.Consumption.City..L.100.km.  Fuel.Consumption.Hwy..L.100.km. 
##                                0                                0 
## Fuel.Consumption.Comb..L.100.km.                             Make 
##                                0                                0 
##                            Model                    Vehicle.Class 
##                                0                                0 
##                        Cylinders                     Transmission 
##                                0                                0 
##                        Fuel.Type 
##                                0

Tidak ada data hilang pada setiap kolom.

1.3 Eksplorasi data

par(mfrow = c(2,6))
hist(data[,1], main="Y - Emisi CO2", cex.main=4)
hist(data[,2], main="X1 - Engine Size.L", cex.main=4)
barplot(table(data[,2]), main="X1 - Engine Size.L", cex.main=4)
hist(data[,3], main="X2 - Fuel Consum.L/km", cex.main=4)
hist(data[,4], main="X3 - Fuel Consum.L/km", cex.main=4)
hist(data[,5], main="X4 - Fuel Consum.L/km", cex.main=4)
barplot(table(data[,6]), main="X5 - Make", cex.main=4)
barplot(table(data[,7]), main="X6 - Model", cex.main=4)
barplot(table(data[,8]), main="X7 - Vehicle.Class", cex.main=4)
barplot(table(data[,9]), main="X8 - Cylinders", cex.main=4)
barplot(table(data[,10]), main="X9 - Transmission", cex.main=4)
barplot(table(data[,11]), main="X10 - Fuel.Type", cex.main=4)

1.3.1 X1 & X8

Dari sini kita bisa melihat bahwa X1 walaupun satuannya liter, namun bukan peubah kontinu. Ini bisa dilihat dari histogramnya yang penuh dnegan celah, mari kita lihat lebih dalam:

unique(data[,2])

## [1] 2 4 6 5 3 7 1 8

Terlihat bahwa jumlah nilai diskritnya terbatas dan sering mewakili kategori yang dibulatkan. Misalnya, mesin dengan ukuran sebenarnya 4.8 liter mungkin dikategorikan sebagai 5 liter.

Selain itu, nilai-nilai ini sering digunakan untuk mengklasifikasikan kendaraan ke dalam kelompok ukuran mesin tertentu, seperti mesin kecil, sedang, dan besar. Oleh karena itu, perlakuan data ini sebagai kategorik lebih tepat dalam beberapa analisis dan visualisasi, mencerminkan penggunaan praktis dalam klasifikasi ukuran mesin.

unique(data[,9])

## [1]  4  6 12  8 10  3  5 16

Begitu pula dengan X8 yakni banyaknya silinder, sehingga kedua peubah ini akan dijadikan sebagai peubah kategorik.

1.3.2 X5 - X9

Saking banyaknya kategori yang ada, bar chart yang dihasilkan nyaris nampak seperti histogram. Mari kita liat:

cat(  "Jumlah baris       =", nrow(data),
    "\nBanyaknya kategori:",
    "\nX5 (Make)          =", length(unique(data[,6])),
    "\nX6 (Model)         =", length(unique(data[,7])),
    "\nX7 (Vehicle.Class) =", length(unique(data[,8])),
    "\nX9 (Transmission)  =", length(unique(data[,10]))
    )

## Jumlah baris       = 7385 
## Banyaknya kategori: 
## X5 (Make)          = 42 
## X6 (Model)         = 2053 
## X7 (Vehicle.Class) = 16 
## X9 (Transmission)  = 27

Terlihat bahwa walaupun bar chart X6 terlihat sangat chaos, namun banyaknya kategori hanya 2000 saja., tidak sama dengan jumlah baris. Ini artinya X6 tidak semuanya unik. Ya sebenarnya bisa langsung dilihat dari tinggi rendahnya bar chart sih. Ini mah biar bisa liat begitu banyaknya kategori model pada X6 saja.

Akan tetapi, random forrest tidak mampu memodelkan lebih dari 53 kategori, sehingga alangkah baiknya jika X6 dikeluarkan dari analisis.

1.3.3 Y (Respon)

Mari kita lihat lebih dalam untuk sebaran data

library(ggplot2)
library(cowplot)

# Histogram
p1 <- ggplot(data, aes(x = `CO2.Emissions.g.km.`)) +
  geom_histogram(binwidth = 20, fill = "#89b49e", color = "black") +
  ylab("Frekuensi") + theme_void()

# Boxplot
p2 <- ggplot(data, aes(y = `CO2.Emissions.g.km.`)) +
  geom_boxplot(fill = "#89b49e", color = "black", outlier.size = 3) +
  xlab("Nilai") + theme_void() + coord_flip()

# Gabung plot
combined_plot <- plot_grid(p1, p2, ncol = 1, align = 'v',
                           rel_heights = c(1, 0.25))

# Menambahkan judul utama
title <- ggdraw() + 
  draw_label("Y - Emisi CO2", fontface = 'bold', x = 0.5, hjust = 0.5, size = 20)

# Menampilkan plot dengan judul utama
plot_grid(title, combined_plot, ncol = 1, rel_heights = c(0.1, 1))

Terlihat bahwa data menjulur ke kanan dan nampak memiliki banyak sekali outlier di sisi kanan.

1.4 Cleaning Data

Ubah tipe data

data <- data[,-7] %>% # Mengeluarkan X6
  mutate_at(vars(2, 9), as.character) # Mengubah menjadi character 

data <- data %>%
  mutate_if(is.character, as.factor) # Mengubah char menjadi faktor
str(data)

## 'data.frame':    7385 obs. of  10 variables:
##  $ CO2.Emissions.g.km.             : int  196 221 136 255 244 230 232 255 267 212 ...
##  $ Engine.Size.L.                  : Factor w/ 8 levels "1","2","3","4",..: 2 2 2 4 4 4 4 4 4 2 ...
##  $ Fuel.Consumption.City..L.100.km.: int  10 11 6 13 12 12 12 13 13 11 ...
##  $ Fuel.Consumption.Hwy..L.100.km. : int  7 8 6 9 9 8 8 9 10 8 ...
##  $ Fuel.Consumption.Comb..L.100.km.: int  9 10 6 11 11 10 10 11 12 9 ...
##  $ Make                            : Factor w/ 42 levels "ACURA","ALFA ROMEO",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Vehicle.Class                   : Factor w/ 16 levels "COMPACT","FULL-SIZE",..: 1 1 1 12 12 3 3 3 3 1 ...
##  $ Cylinders                       : int  4 4 4 6 6 6 6 6 6 4 ...
##  $ Transmission                    : Factor w/ 27 levels "A10","A4","A5",..: 15 26 23 16 16 16 16 16 26 15 ...
##  $ Fuel.Type                       : Factor w/ 5 levels "D","E","N","X",..: 5 5 5 5 5 5 5 5 5 5 ...

Dengan ini data sudah benar dan siap untuk dianalisis.

1.5 Split Data

library('caret')
set.seed(1401211016) # karena partisi data itu random, seperti sampling, maka perlu set.seed agar tidak berubah

# Membuat index partisi data
data.a <- data
colnames(data.a) <- c("Y","X1","X2","X3","X4","X5","X7","X8","X9","X10")
trainIndex <- createDataPartition(data.a$Y, # data Y
                                  p = .8, # 80% data Train, 20% data Test
                                  list = FALSE)
# Membagi data train dan test
TRAIN <- data.a[trainIndex, ] # ambil data train aja
TEST <- data.a[-trainIndex, ] # ambil bukan data train

2 Modeling

K-fold cross validation adalah salah satu teknik validasi untuk mencari tuning parameter terbaik sekaligus mengevaluasi kinerja model. Pada studi kasus ini digunakan 5-fold cross validation. Data dipartisi secara acak ke dalam lima subset data. Secara bergantian masing-masing subset akan dijadikan sebagai data testing, sementara empat subset data lainnya sebagai data training.

library('randomForest')
fitControl <- trainControl(
  method = "cv", # metode K-fold cross validation
  number = 5, # menggunakan 5-fold
  returnResamp = "all")

2.1 Train model

2.1.1 Tuning Parameter dengan `tuneLength`

Opsi tuneLength pada fungsi caret::train akan memilih sejumlah tuning parameter atau kombinasi tuning parameter yang dianggap paling tepat sesuai dengan metode yang dipilih dan data training yang diberikan.

2.1.1.1 Cross Validation

library('randomForest')
rf <- train(Y ~ ., 
            data = TRAIN,
            method = 'ranger',
            tuneLength = 10, # mencoba 10 kombinasi tuning parameter
            importance = "impurity",
            trControl = fitControl, # K-fold cross validation
            verbose = FALSE)
saveRDS(rf, file = "rf1.rds") # Save model di PC
rf

## Random Forest 
## 
## 5908 samples
##    9 predictor
## 
## No pre-processing
## Resampling: Cross-Validated (5 fold) 
## Summary of sample sizes: 4728, 4727, 4725, 4727, 4725 
## Resampling results across tuning parameters:
## 
##   mtry  splitrule   RMSE       Rsquared   MAE      
##    2    variance    26.059115  0.8917964  17.950466
##    2    extratrees  31.196727  0.8708193  23.809074
##   12    variance     6.264561  0.9889481   4.044398
##   12    extratrees   7.282406  0.9855205   4.915968
##   23    variance     5.069730  0.9925003   3.348303
##   23    extratrees   5.368454  0.9916653   3.496834
##   33    variance     4.839064  0.9931297   3.264002
##   33    extratrees   4.988263  0.9927318   3.299288
##   44    variance     4.552820  0.9939448   3.192870
##   44    extratrees   4.877054  0.9930225   3.232910
##   54    variance     4.728557  0.9934049   3.240143
##   54    extratrees   4.805850  0.9932157   3.211387
##   65    variance     4.707321  0.9934568   3.242472
##   65    extratrees   4.736722  0.9933976   3.200456
##   75    variance     4.704185  0.9934591   3.240592
##   75    extratrees   4.727946  0.9934194   3.198672
##   86    variance     4.709383  0.9934453   3.242050
##   86    extratrees   4.712715  0.9934574   3.193748
##   97    variance     4.726564  0.9933934   3.249581
##   97    extratrees   4.694268  0.9935061   3.194368
## 
## Tuning parameter 'min.node.size' was held constant at a value of 5
## RMSE was used to select the optimal model using the smallest value.
## The final values used for the model were mtry = 44, splitrule = variance
##  and min.node.size = 5.

Hasil validasi silang ditampilkan pada plot berikut:

plot(rf, main = "5-Fold Cross Validation Random Forest: tuneLength")

# Selain plot, bisa menggunakan bestTune
rf_best <- rf$bestTune
rf_best

##   mtry splitrule min.node.size
## 9   44  variance             5

Berdasarkan output di atas, tuning parameter terbaik adalah mtry = 44, splitrule = variance dan min.node.size = 5, yang memberikan RMSE = 4.552820, R-squared = 0.9939448 dan MAE = 3.192870.

2.1.1.2 Re-Fit Model Menggunakan Tuning Parameter Terbaik

Jadi setelah mencobakan kombinasi tuning parameter dan di dapatkan tuning, terbaik. Maka selanjutnya kita kan modeling dengan tuning parameter tersebut.

Re-fit model terhadap seluruh data testing dengan menggunakan tuning parameter terbaik yang diperoleh pada tahap sebelumnya:

rf <- train(Y ~ ., 
            data = TRAIN,
            method = 'ranger',
            tuneGrid  = rf_best, 
            importance = "impurity",
            verbose = FALSE)
saveRDS(rf, file = "rf2.rds")
rf

## Random Forest 
## 
## 5908 samples
##    9 predictor
## 
## No pre-processing
## Resampling: Bootstrapped (25 reps) 
## Summary of sample sizes: 5908, 5908, 5908, 5908, 5908, 5908, ... 
## Resampling results:
## 
##   RMSE      Rsquared   MAE     
##   4.817988  0.9931638  3.301318
## 
## Tuning parameter 'mtry' was held constant at a value of 97
## Tuning
##  parameter 'splitrule' was held constant at a value of extratrees
## 
## Tuning parameter 'min.node.size' was held constant at a value of 5

rf_result <- rf$results
rf_result

##   mtry  splitrule min.node.size     RMSE  Rsquared      MAE    RMSESD
## 1   97 extratrees             5 4.817988 0.9931638 3.301318 0.2288119
##     RsquaredSD      MAESD
## 1 0.0006747375 0.06156657

Diperoleh model dengan RMSE = 4.817988, R-Squared = 0.9931638, MAE = 3.301318.

2.1.1.3 Evaluasi Terhadap Data Test

Untuk menguji kinerja model dalam memprediksi data baru, dilakukan evaluasi terhadap data testing:

library(MLmetrics)
eval_test_data <- function(model){
  pred <- predict(model, newdata = TEST)
  mae <- MAE(TEST$Y, pred)
  rmse <- RMSE(TEST$Y, pred)
  R2 <- R2_Score(TEST$Y, pred)
  return(c(RMSE = rmse, R_Squared = R2, MAE = mae))
}

rf_eval <- eval_test_data(rf)
rf_eval

##      RMSE R_Squared       MAE 
##  4.508533  0.994106  3.154478

Diperoleh model RMSE = 4.508533, R-Squared = 0.994106, MAE = 3.154478.

2.1.1.4 Variable Importance

plot(varImp(rf), 
     main = "Random Forest Variable Importance" )

Berdasarkan output di atas, tiga peubah terpenting adalah X4, X2, dan X3.

2.1.2 Tuning Parameter dengan `tuneGrid`

Opsi tuneGrid pada fungsi caret::train memberikan keleluasaan kepada analis untuk menentukan kandidat tuning parameter atau kombinasi tuning parameter yang akan diuji.

2.1.2.1 Cross Validation

# untuk tuneGrid. Jadi kalo tuneLenght itu sangat random, hanya memasukkan angka saja. Sedangkan untuk tuneGrid itu bisa diatur untuk mtru, splitrule dan min.node.size nya
tg <- expand.grid(
  mtry = seq(2, 14, 2), # mencoba kombinasi mtry
  splitrule = c("variance","extratrees"),
  min.node.size = c(5, 10, 20, 30))

rf_tg <- train(Y ~ ., 
            data = TRAIN,
            method = 'ranger',
            tuneGrid = tg,
            ntree = 500, # banyaknya pohon
            max_deep = 100, # maksimal kedalaman
            importance = "impurity",
            trControl = fitControl,
            verbose = FALSE)
saveRDS(rf_tg, file = "rf_tg1.rds")
rf_tg

## Random Forest 
## 
## 5908 samples
##    9 predictor
## 
## No pre-processing
## Resampling: Cross-Validated (5 fold) 
## Summary of sample sizes: 4726, 4727, 4727, 4726, 4726 
## Resampling results across tuning parameters:
## 
##   mtry  splitrule   min.node.size  RMSE       Rsquared   MAE      
##    2    variance     5             25.830164  0.8934379  17.812224
##    2    variance    10             26.073955  0.8923329  18.048965
##    2    variance    20             25.844875  0.8931036  17.842438
##    2    variance    30             26.104912  0.8900580  18.019345
##    2    extratrees   5             31.043809  0.8723268  23.777771
##    2    extratrees  10             30.904549  0.8715876  23.630935
##    2    extratrees  20             30.923607  0.8717697  23.645912
##    2    extratrees  30             31.071931  0.8716321  23.781352
##    4    variance     5             14.711726  0.9515277   9.680727
##    4    variance    10             14.773548  0.9511933   9.697261
##    4    variance    20             14.935770  0.9499941   9.764463
##    4    variance    30             15.080287  0.9489818   9.875269
##    4    extratrees   5             18.165992  0.9331870  13.482297
##    4    extratrees  10             18.285840  0.9322758  13.544688
##    4    extratrees  20             18.531676  0.9305220  13.679850
##    4    extratrees  30             18.368938  0.9315026  13.566974
##    6    variance     5             10.168935  0.9742071   6.641120
##    6    variance    10             10.267738  0.9737277   6.700751
##    6    variance    20             10.581552  0.9719700   6.854151
##    6    variance    30             10.911774  0.9702737   7.014029
##    6    extratrees   5             12.777425  0.9618376   9.343349
##    6    extratrees  10             12.918509  0.9609311   9.427109
##    6    extratrees  20             13.100203  0.9596644   9.472225
##    6    extratrees  30             13.506931  0.9571266   9.761551
##    8    variance     5              7.900230  0.9835605   5.212337
##    8    variance    10              8.096688  0.9827023   5.283437
##    8    variance    20              8.481752  0.9810353   5.497361
##    8    variance    30              8.954909  0.9789076   5.721117
##    8    extratrees   5              9.807478  0.9758535   7.019274
##    8    extratrees  10             10.056662  0.9746185   7.178155
##    8    extratrees  20             10.401438  0.9727993   7.370581
##    8    extratrees  30             10.971568  0.9698289   7.741080
##   10    variance     5              6.705075  0.9876887   4.452409
##   10    variance    10              6.852877  0.9871684   4.537926
##   10    variance    20              7.359380  0.9852425   4.786571
##   10    variance    30              7.827385  0.9833722   5.026027
##   10    extratrees   5              8.130569  0.9825692   5.735172
##   10    extratrees  10              8.380565  0.9815241   5.861804
##   10    extratrees  20              8.858418  0.9795298   6.154391
##   10    extratrees  30              9.428799  0.9768549   6.501314
##   12    variance     5              6.013781  0.9898478   4.037310
##   12    variance    10              6.217857  0.9891914   4.142780
##   12    variance    20              6.680086  0.9875902   4.378622
##   12    variance    30              7.208772  0.9856639   4.635387
##   12    extratrees   5              7.045223  0.9865077   4.879894
##   12    extratrees  10              7.304727  0.9855598   5.040823
##   12    extratrees  20              7.883341  0.9833140   5.403435
##   12    extratrees  30              8.460683  0.9809012   5.772302
##   14    variance     5              5.583719  0.9911311   3.766388
##   14    variance    10              5.799456  0.9904454   3.884107
##   14    variance    20              6.311742  0.9887875   4.146248
##   14    variance    30              6.764581  0.9872051   4.391941
##   14    extratrees   5              6.310069  0.9889568   4.346595
##   14    extratrees  10              6.616041  0.9878927   4.524314
##   14    extratrees  20              7.196900  0.9858214   4.883709
##   14    extratrees  30              7.769542  0.9836235   5.250298
## 
## RMSE was used to select the optimal model using the smallest value.
## The final values used for the model were mtry = 14, splitrule = variance
##  and min.node.size = 5.

plot(rf_tg, main = "5-Fold Cross Validation Random Forest: tuneGrid")

Lihat mana yang RMSE nya paling rendah, itu yang paling baik. Dari sini terlihat bahwa yang terbaik itu adalah variance.

rf_tg_best <- rf_tg$bestTune
rf_tg_best

##    mtry splitrule min.node.size
## 49   14  variance             5

Berdasarkan output di atas, tuning parameter terbaik adalah mtry = 14, splitrule = variance dan min.node.size = 5, yang memberikan CV RMSE = 5.583719, R-squared = 0.9911311, dan MAE = 3.766388.

2.1.2.2 Re-Fit Model Menggunakan Tuning Parameter Terbaik